﻿#ifndef _MULTIGRID_H_
#define _MULTIGRID_H_

#include "Eigen/Eigen"
#include <math.h>

#include "Function.h"
#include "Interpolation.h"
#include "basis.h"

using namespace std;

struct runningLog
{
    double err;
    int count;

    runningLog();
    runningLog(double err, int count = 0);
    runningLog(const Eigen::VectorXd &err, int count = 0);

    inline double updateError(double newErr)
    {
        err = newErr;
        return err;
    }

    inline double updateError(const Eigen::VectorXd &err)
    {
        return updateError(vectorXdInftyNum(err));
    }

    inline const runningLog operator++(int)
    {
        runningLog old(err, count);
        count++;
        return old;
    }

    inline runningLog &operator++()
    {
        count++;
        return *this;
    }

    inline int operator=(const runningLog &l)
    {
        return count;
    }
};

class Multigrid
{
private:
    int degree;
    Function f;
    Dimension dim;
    CycleChoice cycletype;
    Operators::Restriction rtype;
    Operators::Interpolation itype;
    BoundType btype;
    restriction_fun restriction;
    prolongate_fun prolongate;
    Multigrid() = delete;
    void updateDim();

public:
    int size();
    Multigrid(int degree, Function f, Dimension d = One);
    void setDimension(Dimension d);
    void setCycleType(CycleChoice c);
    void setRestrictionType(Operators::Restriction s);
    void setInterpolationType(Operators::Interpolation t);
    void setBoundtype(BoundType b);

    Eigen::SparseMatrix<double> setgrid(int size, Eigen::VectorXd *pValues = nullptr);
    Eigen::VectorXd V_cycle(const Eigen::VectorXd &v, const Eigen::VectorXd &f, int n1, int n2, double w = 2. / 3);
    Eigen::VectorXd FMG(const Eigen::VectorXd &f, int n1, int n2, double w = 2. / 3);
    std::tuple<Eigen::VectorXd, Eigen::VectorXd> Possion_BVP_Multigrid(int nu1, int nu2, int max_iter = 1000, double episilon = 1e-15, double w = 2. / 3, _Out_ runningLog *pLog = nullptr);
};

bool operator==(std::pair<double, double> p1, std::pair<double, double> p2);

Insqaure isinsquare(double x, double y);

inline Insqaure isinsquare(std::pair<double, double> p);

inline double distance(std::pair<double, double> p1, std::pair<double, double> p2);

#endif